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THE BISPECTRUM OF REDSHIFTED 21-CM FLUCTUATIONS FROM THE DARK AGES 
Annalisa Pillepich\ Cristiano Porciani^ and Sabino Matarrese^ 

Abstract 

Brightness-temperature fluctuations in the redshifted 21-cni background from the cosmic dark ages 
are generated by irregularities in the gas-density distribution and can then be used to determine the 
statistical properties of density fluctuations in the early Universe. We first derive the most general 
expansion of brightness-temperature fluctuations up to second order in terms of all the possible sources 
of spatial fluctuations. We then focus on the three-point statistics and compute the angular bispectrum 
' of brightness-temperature fluctuations generated prior to the epoch of hydrogen reionization. For 

I simplicity, we neglect redshift-space distortions. We find that low- frequency radio experiments with 

^— ^ . arcmin angular resolution can easily detect non-Gaussianity produced by non-linear gravity with high 

' signal-to-noise ratio. The bispectrum thus provides a unique test of the gravitational instability 

I scenario for structure formation, and can be used to measure the cosmological parameters. Detecting 

the signature of primordial non-Gaussianity produced during or right after an inflationary period is 
more challenging but still possible. An ideal experiment limited by cosmic variance only and with 
an angular resolution of a few arcsec has the potential to detect primordial non-Gaussianity with a 
non-linearity parameter of /nl 1- Additional sources of error as weak lensing and an imperfect 
foreground subtraction could severely hamper the detection of primordial non-Gaussianity which will 
I benefit from the use of optimal estimators combined with tomographic techniques. 

^ ' Subject headings: cosmology: theory - diffuse radiation - intergalactic medium 
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1. INTRODUCTION 

During the dark ages (the time between recombination and the formation of the first stars), the cosmic microwave 
background (CMB) is coupled to atoms of neutral hydrogen through spin-flip 21-cm transitions. Due to the resonant 
, nature of the interaction, neutral hydrogen at redshift z imprints a signature at a wavelength of 21.12 (1 -I- z) cm in 
■ the CMB. The brightness temperature of the CMB at radio wavelengths thus probes the three-dimensional neutral 
' hydrogen distribution at 30 < z < 100. This accurately traces dark- matter inhomogeneities down to the Jeans 
length (~ 10 comoving pc corresponding to angular separations of ~ 10 ^ arcsec). On smaller scales, the finite 
Q ■ pressure of the gas keeps the baryons uniformly distributed. The expected power-spectrum of fluctuations in the 
?H ' 21-cm background from the dark ages has been recently calculated (Loeb & Zaldarriaga 2004; Bharadwaj & Ali 
2004). In this paper, we focus on the three-point statistics. In particular, we compute the bispectrum of 21-cm 
d maps. This is obtained by expanding the brightness temperature of the cosmic background up to second order in 
the underlying density fluctuations. Related work has been presented by Cooray (2005) who discussed the large-scale 
non-Gaussianity generated during the era of reionization and by Saiyad, Bharadwaj & Pandey (2006) who focussed 
on the post-reionization era (z < 6). 
^ ■ Measuring the bispectrum of fluctuations in the 21-cm background will ultimately allow to quantify the degree 
of non-Gaussianity (NG) of the dark matter density field at high-redshift. Such a non-Gaussianity has two possible 
origins: non-linearity due to the usual Newtonian gravitational instability and non-Gaussianity which is intrinsic in the 
mechanism generating the primordial seeds, i.e. arising during or immediately after inflation, as well as post-Newtonian 
terms inherent in the solution of Einstein's equations to second order (Bartolo, Matarrese & Riotto 2005). The latter 
contributions to non-Gaussianity are usually expressed through a non-linearity parameter /nl which measures the 
strength of quadratic terms in the primordial gravitational potential (e.g. Salopek & Bond 1990, 1991; Gangui et al. 
1994; Verde et al. 2000a). The most stringent limits on /nl are based on the analysis of CMB anisotropics from 
three-year WMAP data (Spergel et al. 2006) and give —54 < /nl < 114 (3(t confidence level). The high-resolution 
observations of CMB anisotropies by the Planck satellite should reduce the 3a detection threshold to |/nl| ~ 10. 
Searches for primordial NG based on the large-scale structure of the Universe where the relevant observables probe 
the dark matter density field (e.g. three-point statistics of the galaxy distribution) are plagued with two problems. 
First, the density field is related to the gravitational potential through the Poisson equation so that its Fourier modes 
weigh the NG contributions with extra terms and only the largest scales may keep memory of primordial NG. 
Moreover, the ratio of Newtonian non-linear terms to intrinsical NG terms roughly scales like (1 -I- z)~^, so that NG 
signatures in the local Universe are easily masked by the effects of gravitational instability. In this sense, the 21-cm 
background at large redhifts appears as an extremely promising dataset to search for primordial NG able to provide 
constraints on /nl which can be complementary and possibly competitive with those based on CMB anisotropies. 

Low-frequency radio observations with high angular resolution are limited by the Earth's ionosphere. This layer 
of the atmosphere is opaque to electromagnetic waves whose frequency lies below the local plasma frequency which 
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corresponds to ~ 15 MHz during day-time and near sunspot maximum and to ~ 10 MHz at night near sunspot 
minimum. These figures are somewhat lower (down to ~ 2 MHz) at preferred sites located near the magnetic poles 
but, in general, the propagation properties of radio waves with frequencies ^' < 30 MHz show extreme variations due 
to ionospheric scintillation (e.g. Kassim et al. 1993). Therefore, it is unlikely that ground-based experiments with 
high angular resolution (< 1 arcmin) will ever measure the properties of the 21-cm background at a redshift z ^ 50. 
Space observatories are needed to overcome these limits. 

The layout of the paper is as follows. In SJ2]we summarize the physics of 21-cm radiation and perform the expansion 
of its brightness temperature to the second perturbative order in all the relevant quantities. Non-Gaussian density 
fields arc introduced in The calculation of the angular bispectrum of redshifted 21-cm fluctuations is presented in 
21 Finally, in SJSJ we estimate the detectability of non-Gaussian features (both primordial and gravity induced) by 
future experiments and conclude. 

We consider a "concordance" cosmological model with a present-day matter-density parameter Jim = 0.27 (of 
which rib — 0.049 in baryons), a cosmological costant contribution of Ha = 0.73 and a Hubble constant of — 
lOO/ikms^^Mpc"^ with h = 0.71. We also assume a Harrison- Zel'dovich primordial power spectrum with spectral 
index n = 1 and a cold dark matter (CDM) transfer function with present-day normalization = 0.9 (where CTg 
denotes the linear rms fluctuation within a sphere with a comoving radius of 8 h^^ Mpc). 

2. THE 21-CM FLUCTUATIONS 

2.1. Kinetic temperature of the IGM 

About 3 X 10^ yr after inflation, atomic nuclei and electrons combined to form neutral atoms, and greatly reduced 
their coupling with photons. Soon after, when the free-streaming CMB radiation cooled below 3000 K and shifted 
into the infrared, the Universe became dark. The cosmic dark ages persisted until the first luminous sources formed 
inside virialized density perturbations [z ~ 30), and reionized the opaque IGM. The exact timing of this event is still 
uncertain, even though observational data bracket it between redshift 20 and 6 (e.g. Loeb & Barkana 2001; Spergel 
et al. 2006). 

The kinetics of recombination and ionization processes in a rapidly expanding Universe leaves a residual degree of 
ionization even during the dark ages. At z ~ 200, the mean hydrogen ionization fraction x ~ 10^** while the fraction 
of singly ionized Helium is much lower (~ 10~®). These figures slightly decrease with time due to the cosmic expansion 
but always remain of the same order of magnitude (e.g. Seager et al. 2000 and references therein). 

Compton scattering of CMB photons off free electrons couples the kinetic temperature of cosmic gas Tgas with the 
temperature of the photon background (Weymann 1965) 

^Tgas 2 Tgas ^Tltot ^ U 

at 3 ntot ot 6 nicC ntot 

where ctt is the Thomson scattering cross section, the electron mass and c the speed of light in vacuum. The total 
number density of particles, ntot, appears on the right-hand side because collisions and Coulomb scatterings hold all 
the different species (electrons, ions and atoms) at the same temperature. Here U represents the radiation energy 
density U = Ui, dv where = A'kJi,/c and — (l/47r) 1,^ dh with 1^ the specific intensity of radiation (energy 
per unit time, surface, frequency and solid angle) and dfi the solid-angle element. In thermal equilibrium, the radiation 
field has a frequency distribution given by the Planck function = B^ITcmb) and the total energy is given by Stefan's 
law U = (tt^ fcg/lS c^) T^]y[g where and h denote the Boltzmann and the Planck constants, respectively. Note 
that the spectrum of the CMB remains close to a blackbody because the heat capacity of radiation is very much larger 
than that of matter (i.e. there are vastly more photons than baryons). For a uniform distribution of baryons with 
physical number density ntot oc (1 -I- z)'^, Compton scattering is efficient down to a redshift z ~ 150 (r2b^^/0.023)^/^ 
(the decouphng era, e.g. Peebles 1993) and keeps Tgas = Tcmb = 2.73 (1 + z) K. Subsequently, the gas expands 
adiabatically, Tgas oc {1 + z^. 

2.2. Spin temperature of neutral hydrogen 

Neutral hydrogen interacts with radiation of wavelength X2i~ 21.12 cm through the resonant transition between the 
two hyperfine levels of the Is state. This corresponds to a frequency 1^21 = 1420.4 MHz and a temperature T* = 0.068 
K. It is convenient to express the ratio between the number density of hydrogen atoms in the excited state, n*,and in 
the ground state, no, in terms of a Boltzmann factor which defines the spin temperature, Tg'. 

^ = ^exp(4), (2) 
no go \ Ts J 

where = 3 and go = 1 are the spin-degeneracy factors of the excited and fundamental levels in the Is state. The 
total density of neutral hydrogen atoms is nni = no + n,. During the dark ages, the time evolution of no is the result 
of the combined effect of atomic collisions and radiative interactions with the photon background. Approximating the 
21-cm transition as infinitely sharp, one gets 



(3) 
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where a{t) = (1 + z)^^ is the cosmic scale factor, Aiq = 2.85 x 10^^^ s^^ is the spontaneous decay rate of the hyperfine 
transition of atomic hydrogen, Coi and Cio are the colhsional excitation and de-excitation rates, while Bqi and Bio 
are the Einstein rate coefficients. At equilibrium, detailed balancing of coUisional and radiative processes gives 

^ = ^exp(-^), ^.^4.^._ = (4) 

L-IO go \ J gas/ -DOl 9* J^Ol 9* 

with Fi, = Anhiy^/c^. When Tgas,IcMB and Tg are much larger than T*, all the Boltzmann and the Bose-Einstein 
factors can be linearized, so as to obtain: 

1 + {l/9o) i (I) ^ (^t " I) + " ^) ' 

where the colhsional de-excitation rate, Cio = (4/3) kio tt-hi, with kiq a function of Tgas tabulated by Zygelman (2005). 
At mean density, the characteristic timescales for reaching an equilibrium state are T*/(Aio Tcmb) — 5.4 x 10^ [(1 + 
z)/51]~^ yr for radiative processes and C^q ~ 1.2 x lO'' (1 -I- 5)~^[(1 -I- z)/51]~^ yr for atomic colhsions (assuming that 
Kio Tgas which is accurate for z ~ 50 where Tgas ^ 50 K). Since these timescales are much shorter than the Hubble 
time, for cosmological studies one is only interested in the steady state solution of equation ([5]) 

TbMB + YcT^ y 

1 + Fc ' Aw Tgas ^ ' 

In other words, atomic collisions drive Tg towards Tgas while 21-cm transitions drive it towards Tcmb- Which process 
dominates depends on the local gas density and collisions become less and less effective with time. For regions at mean 
density, Tg departs from Tgas &i z 100 and approaches Tcmb at z 30. On the other hand, at z 50 collisions are 
efScient only within regions with overdensities such that 1 + 5 > 2.2 [(1 -I- z)/51]~'*. 

This picture is modified by the appearance of the first luminous objects at z < 30. At this epoch, bubbles of 
ionized material start appearing around the sources of light while neutral gas is heated by collisions with fast electrons 
produced by penetrating photons (X-rays and UV radiation) . The detailed evolution of both these phenomena depends 
on uncertain astrophysical details about the formation of massive stars. X-ray binaries and accreting black holes. 
Anyway, the presence of Lya photons provides an additional coupling between Tg and Tgas through the so called 
Wouthuysen-Field effect (Wouthuysen 1952; Field 1958). According to this mechanism, the populations of the singlet 
and triplet Is states can be mixed up by the absorption of Lya photons followed by a radiative decay. The efficiency 
of the effect is determined by the frequency dependence and the intensity of the background radiation field Ji, near 
Lya . This is conveniently expressed in terms of the coupling coefficient Ya (the Lya pumping efficiency) and the color 
temperature Tq, defined as 

^ = _aiogA^ 

where TV^ = J^/F^, gives the photon occupation number at frequency v and Pio = (4/27)Pq, is the indirect de- 
excitation rate of the triplet via absorption of a Lya photon to the n — 2 level with Pq, = 47r J a^, {Ji, /2nhv) dv the 
total Lya scattering rate and an the cross section for Lya scattering (Madau et al. 1997). In steady state, the spin 
temperature of neutral hydrogen is a weighted average of Tgas, Tcmb and Tq, 

„ Tcmb + Tgas + Ya Tq 

i + Fc + r. ■ 

There exists then a critical value of Pq = 27 Aiq Tcmb/4 T, ~ 8.5 x lO^^^ [(1 + z)/ll] s^^ above which Tg ~ Tq. This 
roughly corresponds to '--^ 1 Lya photon per baryon. Note that, in general, background photons emitted as continuum, 
non-ionizing UV radiation and redshifted into the Lya resonance represents the dominant source of coupling for Tg with 
respect to line photons locally emitted via recombinations and colhsional excitations (Madau et al. 1997). Moreover, 
due to the large cross section for resonant scattering and the recoil of the atoms during photon emission, thermalization 
of the radiation is very efficient and in virtually all the physically relevant situations Tq ~ Tgas (Field 1959b). 

2.3. Brightness temperature of 21-cm fluctuations 

The brightness temperature Tb of a radiation field is defined as Ii, = By{Th) and, in the Rayleigh- Jeans limit 
(hv jk^Th ^ 1), Tb = /iy/2 i^^fce- In the absence of energy exchanges between matter and radiation after decoupling, 
the brightness temperature of the CMB would be Tb = Tcmb- However, interactions between CMB photons and atoms 
of neutral hydrogen along the line of sight modify the brightness temperature of the CMB at radio wavelengths. The 
optical depth for the hyperfine transition at A = A21 (1 -f z) (once again approximated as infinitely sharp) is (Field 
1959a) 

- 3c^fiAl0 TlHI 

lefeBJ^IiTg (SWr/c^rphys) 

where the number density of neutral hydrogen, the spin temperature and the radial velocity gradient (in physical units) 
are both evaluated at a comoving distance r(z) along the line of sight. This is obtained in the Sobolev approximation 
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to account for line photons redshifted in or out of the interaction frequency range by the combined action of Hubble 
expansion and peculiar velocities (e.g. Sobolev 1960) In the Rayleigh- Jeans hmit, the radiative transfer equation for 
the rest-frame brightness temperature of a patch of the sky at redshift z gives 

Tb = TcMB e-^^i + (1 - . (10) 

The effect of the 21-cm spin- flip transition can then be measured as an excess or a decrement of the brightness 
temperature of the sky (as observed at the present time) at a wavelength A — A21 (1 -I- z) (Scott & Rees 1990): 

Tb — TcMB Is — TcMB , N 

^21 = — Til ■ (11) 

1+z 1+z ^ ' 

The last equality holds in the optically thin limit which always applies since, neglecting peculiar velocities, 



The spatial dependence of T2i,Ts,Tcmb and T21 is understood in equation (jlip but all these quantities should be 
evaluated at a comoving distance r(z) along the line of sight. Note that intrinsic fluctuations of Tcmb can be fully 
neglected in this context since they are of the order of a few /iK while those of T21 are roughly a thousand times 
larger. Given the thermal history of the gas, four regimes can be distinguished: i) the pre-decoupling era (z > 200) 
when Tg = Tcmb and no fluctuations are present in the brightness temperature of the 21-cm background; ii) the 
dark ages (30 ^ z < 200) when Tg ^ Tcmb and neutral hydrogen absorbs CMB photons; Hi) the reionization epoch 
(10 < z < 30) when Tg ^ Tcmb and neutral hydrogen emits CMB photons; iv) the post-reionization era (z < 10) 
when the IGM is nearly fully ionized over most of its volume and T21 = almost everywhere (with the exception of 
the highest density regions). 

2.4. Fluctuations of the brightness temperature 

The brightness temperature of the CMB at A = A21 (1 -I- z) cm fluctuates with the position on the sky and redshift 
reflecting the variation of a number of underlying quantities. An excess in the density of neutral hydrogen (cither 
determined by an overdensity 5 or by a fluctuation of the ionized fraction S^), 

nni =nH [1-5 (1 + 5,)] (1 + 5) , (13) 

produces an increment of the optical depth and alters Tg (directly, by increasing the collision rate, and indirectly, by 
changing the gas temperature). Radial peculiar velocities, Wpcc, also modify T21 through the term 

^^■^ ^ i/(z) + (l + z)^, (14) 



Sfphys / Or 

where H{z) is the Hubble parameter at redshift z. Finally, during reionization, also fluctuations in the intensity and 
temperature of the Lya background become important. 

For small fluctuations, we perform a perturbative expansion of the interesting quantities by writing di = ^„ S^^"^ /n\ 
where 5|"' = 0[{S[^^y^] and the subscript i indicates a particular random field (e.g. temperature, density, etc.). In 
order to compute the first non-vanishing term of the bispectrum we need to expand T21 up to 2"^^ order in 5, Sx, Wpcc, 
and Sa (i.e. the fluctuation in the pumping efficiency of the Lya background, Y^). Fluctuations of Tgas = Tgas (1 + St) 
can be expressed in terms of the density perturbations through the following functional expansion in real space: 

(5T(r, z) = J d^xi t/i(xi, z) (5(r + xi, z) + J d^xi J d^X2 t/2(xi, X2, z) (5(r + xi, z) (5(r + X2, z) -I- ... (15) 

where the dots denote n-point non-local dependencies with n > 2. Both kernels Gi and Q2 vanish before decoupling 
when Compton scattering keeps Tgas = Tcmb- On the other hand, at late times, when the gas expands adiabatically, all 
the terms in equation (fTSjl are local: Qi{x,z) = (2/3)5d(x) and fJ2(xi,X2,z) = — (2/9) (5d(xi) (5d(x2) where 5d(x) is 
the Dirac delta function. The complete redshift evolution between these asymptotic regimes is determined by equation 
((H), that can be rewritten as 

where, C{t) — {SaT^U Tcmb)/(3 cntot Tgas). Solutions for the kernels Qi and Q2 are derived in Appendix [XI 
where it is shown that tJi(x, z) — gi{z) 5d(x) at all epochs while G2 has a complex spatial dependence. These results 
extend the linear analysis by Bharadwaj & Ali (2004) to second order (our function gi corresponds to their linear 
parameter g). Note that these authors use a similar technique to relate fluctuations in Ts and 5 through equation 
Here, instead, we perturb the steady state solution given in equation ([H]) which allows us to explicit the dependence 
on the underlying variables in analytic form. 

The most general expansion of T21 in terms of all the possible sources of spatial fluctuations up to second order is 
reported in Appendix IB] Here, we focus on the 21-cm signal coming from the dark ages. We therefore neglect fluctua- 
tions in the Lya background and we consider inhomogeneities in the hydrogen ionization fraction whose contribution 
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Fig. 1. — Left: Redshift evolution of the coefficients that couple the brightness-temperature fluctuations to perturbations in the gas 
temperature and density (see the main text and both the Appendices for details). Right: Evolution of the "effective" coupling coefficients 
/li /2 and /s = f"^ introduced in equation l lf 71 1. 

is depressed by a factor of ^ x/{l — xY~ 10~* with respect to that generated by density fluctuations. As a first step, 
we also ignore peculiar velocities, postponing to future work a detailed study of their subdominant effects. In this 
case, the spatial fluctuations of the brightness temperature can be written as 



Ti°)(z) 



Ar2i(r,z)=T2i(r,z) 

= /i(z) (^S^'Hr,z) + ^S<^'\r,z)^+f2{z) S^'\r,zf + 



(17) 



y^^xi J d^X2 tJ2(xi,X2,z)(5(r-|-xi,z)^(r + X2,2;) 



where 



rp{0) _ r 
J 21 - Jo 



, T'cmb\ _x 
1 — j{l-x) 



s 



and, since H{z) ~ Hq nU^ (1 + z)^/^ for z > 1, 
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69.05 mK 
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(18) 



(19) 



\ 16A:bi^21 , 

The functions /i, /2 and /a are different combinations of the coefhcients that couple the brightness temperature with 
the local baryonic density (superscript b) and the gas temperature (superscript T): 

fi^fi+gifl; (20) 



/2=/2 



bb 



9i fr+9i f: 



bT. 



(21) 

fs^fl- (22) 
Analytic expressions for these coupling constants are given in Appendix IB] and their behaviour as a function of redshift 
is shown in Figure [TJ An overdensity of order unity gives rise to a fluctuation in T21 of a few tens of mK and the 
expected signal is maximum at 2; ~ 50. Note that /2 and /a change sign with redshift and multifrequency experiments 
can be used to isolate the contributions of the different terms in equation pT|) . 

3. NON-GAUSSIAN DENSITY FIELDS 
3.1. The bispectrum 

From now on it will be convenient to work in Fourier space. We define the Fourier transform of the random field 
X(r) as X(k) = J cfir X{r) exp (— ik • r), where r = rh. For a stationary random field (statistically invariant under 
rotations and translations), the power spectrum P(fc) is given by the ensemble average 

(JX(ki) dX{k2)) = (2^)3 P(fci) fc(ki + k2) . (23) 

Similarly, the bispectrum i?(ki,k2,k3) is defined as 

(,5l(ki) SXik2) SXiks)) = (27r)3 B(ki, ks, kg) ^(ki + ka + ka) . (24) 

For a Gaussian random field, B = 0. 
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3.2. Non-Gaussianity from gravitational instability 

For the range of scales and redshifts considered here, the baryons trace the dark matter, so that we can speak of an 
unique overdensity S = 5'^'' + ^fJ*-^' +0(15^^^), where the Fourier-space expressions for the first and second-order terms 



read, respectively 
and 



jW(k,z) = Z?H.(^)^«(k) (25) 



l^(2)(k,z) = :|t^ J d^q, J d^q2 Su{^l + ^2-^)lC(^l,^2,z)~s(^\^^)~s^„'\^2) , (26) 

where _D_|_(z) is the hnear grow factor of density fluctuations and (5q^^ is the matter overdensity linearly extrapolated 
till the present epoch (i.e. at z = when D+ = 1). For the density field originated via gravitational instability of 
Gaussian initial conditions, the kernel K. assumes the form (e.g. Fry 1984) 

/CL(kl, k2) ^l + l C0S2 0,, + Uh+'ll) cos (27) 

7 7 2 \ fc2 fci / 

with cos = ki • k2 and k — k/fc. This holds at high-redshift (where the mass density parameter approaches unity). 

3.3. Primordial non-Gaussianity 

For a large class of models for the generation of the initial seeds for structure formation, including standard single- 
field and multi-field inflation, the curvaton and the inhomogeneous reheating scenarios, the level of primordial non- 
Gaussianity can be modeled through a quadratic term in Bardeen's gauge-invariant potential <I>, ^ namely 

<f = <i>L + /NL(<i>^ -($!)) , (28) 
where $l is a Gaussian random field and the specific value of the dimensionless non-linearity parameter /nl depends 
on the assumed scenario (e.g. Bartolo et al. 2004). The kernel 

^(qi,q2,2) = /CL(qi,q2) + /CNL(qi,q2, 2) (29) 

with 

J^(|qi + q2|) n , 1 , J_ 

nqi)T{q2) Wl qy qiq2 
accounts for both the non- linear effects of gravitational instability (/Cl) and the evolved non-Gaussian contribution 
(/Cnl). Here E{z) = Hq flom/ (4c^ D+{z)), and T{k) is the linear transfer function for CDM. 

It is worth stressing that equation even though commonly used, is not generally valid: detailed second-order 
calculations of the evolution of perturbations from the inflationary period to the present time show that the quadratic, 
non-Gaussian contribution to the gravitational potential should be represented as a convolution with a kernel /nl(x, y) 
rather than a product (in full analogy with the quadratic term in equation (|15p ). Indeed, for scales which entered the 
Hubble radius during matter dominance, the calculation of /nl can be performed analytically in the so-called Poisson 
gauge (Bartolo, Matarrese & Riotto 2005) and leads to the Fourier space expression: 

...... 5 77(ki.k2)-102(fc? + /ci) , 

/NL(ki,k2,2:j = - (flni-l) 



/CNL(qi,q2,^) = 6hLE{z) ■! (- + - + — COS f?i2 ) (30) 



3 ' ' 42/c2 

(fc2+fc|)(ki.k2) /4,,_ 



)L2L2jL2 V 7 

1 2 \ 

6 kfk^ 6 (ki • k2)2 



kj + k 



0{E{zf/k^) , (31) 



7 /fc4 7 yfc4 

where the parameter Oni— 1 quantifies the amount of non-Gaussianity produced during (or immediately after) inflation. 
For the simplest case of single-field slow-roll inflation, |ani — 1| <C 1, being of the order of the slow-roll parameters 
(Gangui et al. 1994; Acquaviva et al. 2003; Maldacena 2003), while alternative models, like the curvaton or the 
inhomogeneous reheating scenarios, may easily accomodate much larger values of a^i — 1 (see Bartolo et al. 2004 and 
references therein), such that all momentum-dependent contributions to /nl can be neglected and the non-linearity 
parameter can be effectively approximated with a constant, /nl ~ f (o-ni ~ l)i as assumed in the present analysis (see 
equation ([28]) ). 

From the expressions above, one can easily obtain the (tree-level) bispectrum of density perturbations evaluated at 
three different redshifts as 

B(ki,k2,k3,Zi,Z2,Z3) = 2 /C(ki,k2,zi)7?+(zi)2-D+(22)i^+(23)Po(fci)Po(fc2) +cycl. , (32) 
where Po{k) is the linear power-spectrum of density fluctuation extrapolated till the present epoch. The above 
expressions tell us that the optimal strategy to constrain non-Gaussianity in the large-scale structure of the Universe 
is to use large-scale and/or high-redshifts datasets (e.g. Verde, Heavens & Matarrese 2000b; Verde et al. 2001; 
Scoccimarro, Sefusatti & Zaldarriaga 2004). From this point of view, the 21-cm fluctuations potentially represent a 
very powerful tool to detect non-Gaussianity which could complement the limits obtained from the study of CMB 
anisotropics (e.g. Spergel et al. 2006). 

^ On scales much smaller than the Hubble radius, Bardeen's gauge-invariant potential would reduce to minus the usual peculiar gravi- 
tational potential. 
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4. THE ANGULAR BISPECTRUM OF 21-CM FLUCTUATIONS 

4.1. Definition 

The observed brightness temperature is a convolution of T21 with a frequency-dependent instrumental response 
W^(r) (from now on the window function): 

T2i''{n) = Jdr W,{r) nirn). (33) 

The window function peaks at the central frequency of the detector, v, and is different from zero within a characteristic 
frequency range corresponding to the bandwidth of the instrument. For analytic convenience we prefer here to express 
Wi, in terms of the comoving distance r evaluated at z = V2i/v — 1. 
The fluctuations of the observed ^ brightness temperature can be expanded in spherical harmonics: 



Ar2i(n)=^ ^ a,,nYt^{h), a,^ ^ / d^n Ar2i(n) y/„(n) . (34) 

1=0 m=-t 

The angular bispectrum is then given by the ensemble average 

which for an isotropic random field can be re-written in terms of the angle-averaged bispectrum Bi-^i^i^ as 

Rmi^2m3 h h 4 \ ^ (3g) 

«i«2*:3 I mi 1712 ™3 J 1 ^ ' 

Here the matrix denotes the Wigner-3j symbol so that B™^^^'"^ does not vanish only when: 

1. \£j — £k\ < < £j + for all permutations of indices i,j, k ; 

2. £1 + ^2 + 4 = even ; 

3. TOi + TO2 + ™3 — . 

4.2. Derivation 
The angular bispectrum of the 21-cm background is then ^ 

= J dhi dh2 dhs J dn dr2 dr^ W,{n)WAr2)W,ir3) (37) 

(AT2i(ri) Ar2i(r2) Ar2i(r3)> r;,„,(ni) r,;„,,(n2) F;3„3(n3) , 
where the three-point correlation function is 

^-^^ ^ B2l(ki,k2,k3,Zl, 22,23) (38) 

(5D(ki + k2 + k3)e'('^i""^+''^''^+''^""^) 
with i32i(ki, k2, ka, zi, Z2, 23) the three-dimensional bispectrum of AT21. 

4.2.1. Bispectrum from gravitational instability 

We first account for the non-Gaussianity generated by gravitational instability from Gaussian initial conditions. In 
this case, the integral over ks in equation ([38| can be easily computed exploiting the Dirac delta function. We then 
expand the exponential term in the same equation into spherical harmonics and perform all the angular integrations. 
On scales where the baryon distribution traces the density fluctuations of the dark matter, 

B2i(ki, ks, ka, zi, Z2, Z3) = 2 JC2i{k^,k2,zi)D+{zifD+{z2)D+{z3)fi{zi)fi{z2)fi{z3)Po{ki)Po{k2) + cycl. . (39) 

where the kernel /C21 can be written as a power series of cos0i2, 

2/C2i(ki,k2,z) = Ao + Ai cos6'i2-KA2 cos^ 6*12 (40) 

with 

, 10_^/3 , /^_^/3 \fkl^k2\ 4 /3 

Ao^- + jg.. + j^, + ^ A2 = - + -g2.. (41) 

The factors (72a, 52b, 52c are functions of time that parameterize the evolution of the kernel Q2 and thus take into account 
gas temperature perturbations up to second order in S as indicated in Appendix |^ After long but straightforward 

* In order to simplify the notation, from now on we drop the superscript "obs" to indicate the observed brightness temperature. 

^ We just consider the case of a single window function. The generalization of our results to multi-frequency analyses is straightforward. 
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calculations (see Verde et al. 2000b for a similar case), the expression for the angular bispectrum of the 21-cm 
background can be written as a sum of cyclic terms 

Beiezis = Bi^i^ + Bi^i.^ + Bi^e^ (42) 

with 



„ 16 (2^i + l)(2£2 + l)(24 + l) 

^ V (4^) 

ae,{k,)ai,{k2)zPo{k,)P„{k2) ^ z(^«+^^)(-l)^(24 + 1)(2^7 + 1) 



J dki J dk2 kl kl (43) 



(i2hi\ \ili2e3 



oooy\^oooy\^oooy \^74 ^ 

drl?+ (r)2 /i (r) M^, (r) j,, (fcir) j,, (fcar) (fci , , r) , 

where {. . . } indicates a Wigner-6j symbol. Here ai{k) = J drji{kr) D+{r) fi{r) Wi,{r), with ji{x) the spherical Bessel 
function of the first kind. Note that only terms with £q — £1 — £, ..,£1 + i and £y = £2 — £, --,£2 + ^ contribute to the 
sum. The functions Bi{ki, k2, r) are generated by the terms in equation Since this power series does not contain 
exponents larger than 2, only terms with £ = 0, 1, 2 matter. We thus obtain: 

50=^0 + ^^2 = ^ + ^ (^52a + I g2b) + J , (44) 

^2 = 1^2 = 7^ + ^^ .92b. (46) 



3 - 21 3/: 



For narrow window functions, the expression in equation (|43p can be further simplified. In particular, the linear growth 
factors, the coupling constants fi and the time functions g2i introduced in Appendix [Xl can be evaluated at the epoch 
corresponding to the central frequency of the detector. 

4.2.2. Bispectrum from primordial non-Gaussianity 

We follow a different procedure to integrate equations ((37)) and p8|) with the kernel /CNL(qi, q2, z) given in equation 
(j30p . In this case, it is convenient to expand the Dirac delta function in spherical harmonics and perform the angular 
integrations (e.g. Komatsu & Spergel 2001). The result is 



R 16 f / (2^1 + 1)(2^2 + 1)(2£3 + 1) ( h i2 X I , 2 ( .0 ( ( n U7\ 



where 



I3,{r)= j dk^^at{k)j,{kr) and T ,{r) = Q j dkk^ T{k) e,{k) j,{kr) , (48) 

with e,(fc) = / dr' D+{r'Y [E{r') h{r') + g2A{r')h{r')/Q\ W,{r') j,{kr'). Note that, at z ~ 50, 6Efi/g2df3 ^ 20 and 
the contribution of the non-local kernel Q2 to ce is rather small (see also Figure |4|). 

4.3. Results 

We perform a numerical integration of equations (j43p and (|47l) . The resulting angular bispectra are shown in Figure 
[2] for different sets of the parameters £±,£2 and £3 (bispectrum configurations). From top to bottom, we present 
equilateral (Bm), quasi-equilateral (Bn^ie-i) and "squeezed" {Bie2) configurations. In all cases, we consider an 
ideal instrument with unlimited angular resolution and characterized by a Gaussian window function with a frequency 
(Itr) bandwidth of 0.1 MHz centred around z ~ 50. The left panels compare the different contributions due to 
non- linear gravity - equation (|43p - and to primordial non-Gaussianity - equation (|47p - with /nl — 1- The total 
non-Gaussian signal in the various configurations and for different values of /nl is shown in the right panels of FigureO 
Note that the signal due to gravitational instability dominates at all angular scales. For £ < 200, the ratio between the 
two contributions is 50. This difference becomes much more severe at smaller angular scales. This happens because 
any primordial non-Gaussianity in the gravitational potential is going to be shifted to larger scales when observed in 
the overdensity field as a direct consequence of the Poisson equation. 
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10' 10= 10= 10^ 101 10= 10^ 

1 1 



Fig. 2. — The bispectrum of 21-cm anisotropies -B^i.^a.^a measured by an ideal experiment with 0.1 MHz bandwidth centered 
around z = 50. From top to bottom, we consider equilateral (£i =1-1 = ^3), quasi-equilateral (^2 = + 1,^2 = l\ — 1) and 
"squeezed" (^2 = ^1,^3 = 2) configurations, respectively. In the left panels, the contributions of non- linear gravity (solid) and 
of primordial non-Gaussianity (with /nl = 1; dashed) are compared. The right panels show the total bispectrum for different 
values of the non-linearity parameter. Note that in the equilateral and quasi-equilateral configurations the bispectrum vanishes 
when I is odd (not plotted), is positive for i = 2j with j an odd integer and negative for i = 2j with j an even integer. 



10 



The bispectrum of 21-cm fluctuations 



V 




10- 



101 



gravitational 
primordial 

I I I I I 



102 



103 
1 



10* 



105 



V 

m 




10- 



101 



gravitational 
primordial 



102 



103 

1 



104 



105 



Fig. 3. — Cumulative signal-to-noise ratio for tlie measurement of the bispectrum of 21-cm anisotropies using modes up to a maximum 
value of I. Solid lines refer to non-Gaussianity generated by gravity while dashed lines indicate the primordial signal with /nl = 1- These 
quantities have been computed using the signal-to-noise ratio per mode of the quasi-equilateral configuration (left) and of the squeezed 
configuration (right). 

5. DISCUSSION AND CONCLUSIONS 

Using second order perturbation theory, we have computed analytically and evaluated numerically the bispectrum 
of redshifted 21-cm fluctuations down to a few arcsecond angular scales ijl ~ 10^). Is the expected signal detectable 
by future radio experiments? And can the data be used to test the gravitational instability scenario and quantify 
primordial non-Gaussianity? Suppose we fit the observed bispectrum with a theoretical model by minimizing 

the objective function 

2 



X 



E 



(49) 



where i distinguishes the different contributions of non-linear gravity and primordial non-Gaussianity and Pi is a 
normahzation constant (i.e., for a given set of cosmological parameters. Pi = (crg/O-Q)^ and P2 = /nl (o's/O.O)''). If 
the amount of non-Gaussianity is small, the cosmic variance of the bispectrum is given by the disconnected part of 
the six-point function of agm (Luo 1994). The variance of Bg^t^t^ is then calculated as (e.g. Spergel & Goldberg 1999) 



^ili2i3 ~ \Bllt2t3) 



where A^^^^^s is equal to 1, 2 or 6 when all the ^-indices are different, two of them are the same or all of them are 
identical, respectively. Here Cg denotes the total angular power spectrum of the 21-cm background including the 
contribution of the detector noise. The corresponding Fisher matrix Fij is 



E 

ii<i2<e.3 



r)(«) r)(j) 



a 



(51) 



11^2^3 



where i = 1 refers to parameter Pi and i 
z-component is given by 



2 refers to parameter P2. The signal-to-noise ratio {S/N)i for the 
1 



(52) 



An order of magnitude estimation of this quantity as a function of the angular resolution £ can be obtained neglecting 



the covariance of the different components, thus writing {S/N)i ~ ^/Fu. Since the calculation of Bg^g^g^ is expensive at 
high we use the configurations in Figure[2]to estimate the total signal-to-noise ratio. Our results are shown in Figure 
[21 an experiment with arcmin-scale resolution at z = 50 would clearly detect the bispectrum induced by non-linear 
gravity with {S/N) ^ 100. On the other hand, for the signal from primordial non-Gaussianity, (S/N) ~ 0.1 /nl- 
Increasing the angular resolution by a factor of 10 would give (S/N) ~ 0.3 /nl- We consider here an ideal, full-sky 
experiment where the measurements are only limited by cosmic variance (i.e. with no detector noise) and where 
foreground signals can be perfectly subtracted. If only a fraction /sky is surveyed, our signal-to-noise estimates should 



be depressed by a factor / 



1/2 



sky ■ 
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A tomography of the neutral hydrogen distribution within the redshift range 30 < 2: < 100 in slices of bandwidth 
^0.1 MHz (la) would provide nearly 150 21-cm maps. Since the correlations between different redshift slices will be 
negligible for high multiples {£ > r/Ar ~ 1500), this would increase the signal-to-noise ratios by a factor of ~ 10. 

Our results show that studies of the 21-cm bispectrum have the potential to test the gravitational instability scenario 
and the origin of primordial fluctuations. Measurements of three-point statistics will also provide information on the 
cosmological parameters through the shape dependence of Bi-^i^i.^ and the redshift evolution of the growth factors. 

Low-frequency radio observations with high-angular resolution (a few arcsec, corresponding to £ ^ 10^) can detect 
primordial non-Gaussianity with /nl ~ 1- One should however take into account that such high values of ^ correspond 
to perturbation modes which have crossed the Hubble radius during radiation dominance; for such modes and for 
values of /nl ~ 1, a full second-order evaluation of the matter-transfer function would be needed to account for 
e.g. the Meszaros effect beyond the linear approximation, following Bartolo, Matarrese & Riotto (2007). On the 
scales of interest, weak lensing by the large-scale structure of the Universe typically modifies the angular power 
spectrum of 21-cm fluctuations by only ~ 1% (Mandel & Zaldarriaga 2006). However, it also produces a non- vanishing 
disconnected four-point function which represents an additional source of noise for the bispectrum. Potentially this 
could severely hamper the measurement of /nl and a detailed study of the effect is required. Major improvements 
in the measurement of /nl could be achieved by constructing optimal estimators that weigh more those triangular 
configurations that correspond to particularly high signal-to-noise ratios. We will address these issues in future work. 

We thank Nicola Bartolo, Michele Liguori, Simon J. Lilly and Padclis Papadopoulos for helpful discussions. AP 
thanks the organizers of the workshop on non-Gaussianity held on July 2006 in Trieste where the main results of this 
paper have been presented. 
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APPENDIX 

A. PERTURBATIONS IN THE GAS TEMPERATURE AS A FUNCTION OF DENSITY 

FLUCTUATIONS 

Equation [16] shows that temperature and density perturbations in the cosmic distributions of baryons are related 
to each other. In this Appendix we want to express St in terms of S up to second order in perturbation theory. Our 
starting point is the functional expansion given in equation ()15p that we Fourier transform, to obtain: 

~6TiKt) = gi(k,t)S{k,t)+ Jd^'qi y"d3q2 g~2(qi,q2,0fc(qi+q2-k)J(qi,i)^(q2,i) + ... . (Al) 

Performing a perturbative expansion up to second order gives: 

5W(k,i)=e~i(k,t)5(i)(k,i) , (A2) 

i4')(k,<) = igi(k,i)^(2)(k,t) + J d3qg~2(q,k-q,i)^(l)(q,i)^(i)(k-q,i). (A3) 

Thus, Taylor expanding equation (|16p and comparing terms of first perturbative order, one obtains: 

where we factorized the growth of linear overdensities as in equation (j25p . Since this differential equation does not 
depend on k explicitly, it follows that Qi is only a function of time, i.e. t/i(k, i) — gi{t) and tJi(x, t) — gi{t) S-d{'x). 
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Fig . 4. — Redshift evolution and Fourier-space dependence of the kernels Qi and Q2 relating temperature and density perturbations (see 
i|2.4l and Appendix|A]for further details). Top left: redshift evolution of the first-order function gi{z). Top right: contour levels of the kernel 
62(qi,q2,z) at z = 50 with /nl = 1 and qi || q2; consecutive levels increase by a factor of 10 starting from 0.1 at the centre. Bottom 
left: redshift evoluti on o f the functions g2a,92b and 92c introduced in equation I IA6II . Bottom right: redshift evolution of the function g2d 
defined in equation I IA6I I. 



The corresponding function gi{z) can be obtained performing a simple numerical integration and is plotted in Figure 



m 



The evolution equation for second-order fluctuations reads 

dD+it) 



dt 



dgi 
dt 



dgi 

' dt 



.91 CW , 



(A5) 



^(qi,q2,t) 



where we used equations (pS)) . ([77|) . and (|30p to write (5*^^^ in terms of the kernel /C(qi,q2,i)- The solution of 
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this differential equation can be written as follows: 

^2(qi,q2,t)=.92a(<) +92h{t) COS^ 6*12 +52c(i) 



/nl 92d{t) 



r(|qi+q2|) / 1 



-- 








1 


2 



cos 6/12 
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(A6) 



T{qi)T{q2) \qf ' g| 9192 
where cos0i2 = qi •q2- The time-dependent terms (72a, ff2b, 52c and g2d are plotted in Figure|4]as a function of redshift. 

B. SECOND ORDER EXPANSION OF T21 

We present the most general expansion of T21 in terms of all the possible sources of spatial fluctuations up to second 
order. The brightness temperature of the 21-cm background can be written as: 



■ cost/12 



T21 = T^f + 6n{> + - sni 



1 



(2) 



where 



(0) 



/o 1 



(2) 



2 



dr 
1 .„at;(2) 



(Bl) 

(B2) 
(B3) 



if"A(2) 



(B4) 



(9r 



+ /2 



e(l)^ .-x(l)5(l) , /-va^rC 



dr 



dr 



Here 5^ indicates the n-th order fluctuation in a given random field denoted by the subscript i: b stands for the 
baryon density, T for the gas temperature, x for the hydrogen ionization fraction, v for the radial peculiar velocity, 
and a for the Lya pumping efficiency. 

To simplify the notation, we denote zero-th order (i.e. unperturbed), time-dependent quantities with overlined 
symbols and we indicate particular combinations of these zero-th order quantities with a few special characters: 



_ 4 Kio(rgas) 
•J ^10 J gas 



yc = cii + zf{i-x) , 



y 



yc , 



AT = T, 



CMB 



Lgas , 



S — — [TcmB + Ya Tgas + 3^c'7g: 



(B5) 
(B6) 



We also introduce the symbols 771 and r/2 to parameterize the dependence of the collisional coupling rate kiq (defined 
in Section [2T2|) up to second order on Tgas: 



(B7) 



Adopting this notation, the redshift-dependent weights which couple the brightness temperature to the various sources 
of fiuctuations are: 



3c^/l4io\ nm f 1 

16fcBJ^2l/ ^^0 \f^n: 



1/2 



(1+^) 



1/2 



/o 
.fo 



1-^1(1--) 3^T^^ Ar(i 



TcMB C (l + z)3 
3^5 



1 - 2 



AT 
3^' 



TcMB 
3^ 



(>"« Tgas)(l-X) 



(B8) 
(B9) 
(BIO) 
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